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Abstract 

In this paper we classify the allowed order parameter symmetries in multi- 
layer cuprates and discuss their physical consequences, with emphasis on 
Josephson tunneling and impurity scattering. Our solutions to the gap equa- 
tion are based on highly non-specific forms for the inter- and intra-plane pair- 
ing interactions in order to arrive at the most general conclusions. Within 
this framework, the bi-layer (A^ = 2) case is discussed in detail with refer- 
ence to YBCO and BSCCO and the related Landau-Ginzburg free energy 
functional. Particular attention is paid to the role of small orthorhombic 
distortions as would derive from the chains in YBCO and from superlattice 
effects in BSCCO, which give rise to a rich and complex behavior of the multi- 
layer order parameter. This order parameter has N components associated 
with each of the bands or layers. Moreover, these components have spe- 
cific phase relationships. In the orthorhombic bi-layer case the (s, —s) state 
is of special interest, since for a wide range of phase space, this state exhibits 
vr phase shifts in corner Josephson junction experiments. In addition, its 
transition temperature is found to be insensitive to non-magnetic inter-plane 
disorder, as would be present at the rare earth site in YBCO, for example. Of 
particular interest, also, are the role of van Hove singularities which are seen 
to stabilize states with dnc2_y2-\\ke symmetry, (as well as nodeless s-states) 
and to elongate the gap functions along the four van Hove points, thereby 
leading to a substantial region of gaplessness. We find that dri.2 _y2-\Ske states 
are general solutions for repulsive interactions; they possess the fewest number 
of nodes and therefore the highest transition temperatures. In this way, they 
should not be specifically associated with a spin fluctuation driven pairing 
mechanism. 

PACS numbers: 74.72.-h, 74.20.Mn, 74.50. +r, 74.62.-c 
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I. INTRODUCTION 

The question of the order parameter symmetry in layered cuprate superconductors is 
one of the most important issues currently under debate. Despite strong evidence for 
a dx2_y2 state,0 at least in one particular cuprate (YBCO), there are still experimental 
inconsistencies. ili This matter is complicated further by the complexities of YBCO associ- 
ated with the double copper-oxide plane structure. The same bi-layer unit is shared with 
the bismuth 2212 compounds (BSCCO), where recent photoemission data,QS seem also to 
favor a d^^-yi symmetry. In the former compound,! if not in the latter,0 there is evidence 
that the double layer unit leads to two bands each of which crosses the Fermi surface and 
each of which presumably has a distinct superconducting order parameter. It is not a priori 
clear whether these two gaps have the same or opposite phase, nor whether they are predom- 
inantly associated with the same or different irreducible representations of the tetragonal 
symmetry group. It is, therefore, essential to provide a systematic classification of the or- 
der parameter in these bi-layer (and more general A^-layer) systems before meaningful and 
unambiguous inferences can be deduced from the experimental data. 

It is the purpose of this paper to expand upon earlier work by ourselves! and other 
groupsiHii by classifying and establishing physical consequences of the various order pa- 
rameter symmetries in iV-layer svstems. Of particular interest are those parameter sets 
which lead to vr junction behaviorcjil3 in an a, 6-plane Josephson configuration. It is widely 
assumed that the observation of these tt phase shifts provides the strongest evidence yet 
for the dx2-y2 state in YBCO. By contrast, we find that in the more general orthorhombic 
bi-layer case, in part because of van Hove effects, tt phase shifts are fairly widespread and 
not uniquely associated with the c?^2_y2 state. Moreover, the proximity to the van Hove 
singularities leads to a stabilization of the d^2_y2 symmetry for a variety of different models 
for the pairing interaction, beyond the simple spin fluctuation model.lla It also is associated 
with a considerable distortion of the gap function away from the ideal representation of the 
dx2_y2 state. In order to understand the physical consequences of various order parame- 
ter sets we investigate the role of intra- and inter-layer impurities, include a discussion of 
Josephson coupling, and present a more general analysis of our results in the context of both 
a Landau-Ginzburg theory and the solution of the arbitrary A^-layer problem. 

A number of authors have considered the possibility of electron tunneling and pair in- 
teractions of electrons on different planes in the context of a bilayer structure. Bulut and 
Scalapino0 used numerical solutions of a strong coupling model to show that when both 
interlayer interactions and hopping were included in a bilayer model two competing states 
arose, one with dx2_y2 solutions in phase on both sub-bands of the Fermi-surface and the 
other with s-states of opposite sign. Similar observations were made by Liechtenstein and 
co-workersi who argued that the more probable situation for spin fluctuation induced su- 
perconductivity, corresponded to a pair of out-of-phase s states. It should be noted that 
throughout this paper we will use the generic notation "s" ( and "d" ) as applying to gap 
functions which have the same (or different) signs under a rotation of the wave-vector by 
7r/2. It is important to be particularly clear on our notational convention since we emphacize 
the role of orthorhombicity in our work. Among other things, this orthorhombicity leads to 
what we refer to as "s-d" or "mixed states", (although the overall sign under a 7r/2 rota- 
tion remains either + or — , depending on which component is dominant). The possibility 
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of mixed s-(i-states has been discussed within Landau-Ginzburg theory, in a mono-layeJli 
as well as a bi-layer context .0 In the latter, a spontaneous s-d mixing was considered, 
using microscopic arguments based on an RVB decoupling scheme. While most bi-layer 
studies considered tetragonal symmetry, earlier work by our own group,! using the more 
conventional BCS pairing formulation, investigated the role of a small amount of a, 6-axis 
anisotropy in the bi-layer problem. This small orthorhombicity appears to be amplified by 
van Hove effects. As a consequence, the solution corresponding to a pair of out of phase 
s-wave states frequently has the interesting and important physical consequence of leading 
to TT phase shifts in a corner Josephson tunneling experiment. tdiill It should be stressed that 
orthorhombicity is assumed to enter via the chains in YBCO. There are alternate scenarios 
for this material which explicitly build in the chain bands&^ and their contribution to the 
superconducting gap. Here we ignore these explicit effects, in large part because recent tun- 
neling studies indicate that the order parameter behavior is largely unaffected by reducing 
the oxygen stoichiometry to the limit where the chains are strongly fragmentedli^ (and thus 
presumably irrelevant to the superconductivity). 

The more general A^-layer problem has been studied by Bulaevski and Zyskinlll as well 
as Klemm and Liu.El In these calculations, it was also assumed that the c-axis consisted of 
a coherently coupled stack of bi-layer structural units. Here we presume that the consid- 
erable evidence for the absence of c-axis coherenceilii requires a different starting point. 
While, for simplicity, we focus on an isolated A^-layer complex, it is clear that higher order 
effects associated with incoherent coupling between unit cells must ultimately be included. 
Such incoherent coupling can be introduced following, for example. Ref. An additional 
complexity was raised in even earlier studies by Efetov and LarkinE3 who pointed out that 
a complete treatment of interlayer electron pairing required the consideration of a triplet 
pairing state odd under exchange of layer indices. This case was further investigated by 
Kettemann and EfetovEl as well as Klemm and Liu0 who argued that the mixing of triplet 
and singlet states would generally lead to a second transition below T^. Thus far there is no 
evidence for this second transition, nor is there much support for the triplet pairing state. 

There are important issues of controversy implicit in the work discussed here. There 
appears to be little doubt that the components of the bi-layer couple magnetically. However, 
it has not been persuasively demonstrated that there is a coherent electronic coupling t± 
within the unit cell.S Neutron experiments^^ on YBCO provide evidence for a c-axis 
modulation associated with the spacing of the planes within the bi-layer. While this can 
be explained by magnetic coupling alone,iiii earlier photoemission studies on this cuprate 
have reported the observation of two copper oxide plane bands.! The situation in BSCCO 
is even less certain with two photoemission groups reaching opposite conclusions.!'! For the 
purposes of the present paper it will be assumed that there is coherent coupling between 
the layers, though not between the unit cells. This assumption is based in large part on 
the demonstrated intra-bi-layer magnetic interactions which suggest a moderate degree of 
communication between the layers. Moreover, on this basis it may be presumed that there 
are direct or indirect electronic interactions within the bi-layer complex which must be 
included in any complete theory of the superconductivity. 

Similarly controversial is the origin of a possible d^^-^i state. While this state is con- 
sistent with spin fluctuation mediated superconductivity,0il in this paper we demonstrate 
that it appears more generally as a solution to the gap equation in the presence of repulsive 



3 



REVTeX 3.0 



Levin Group Preprint, 1995 



interactions. Providing this repulsion has some wave-vector dependence (near, although not 
necessarily pinned at the anti-ferromagnetic position (7r,7r)), we find the d^^i-y-i state to be 
stable. Thus, experimental observation of this symmetry appears to provide more support 
for superconductivity arising from repulsive interactions than for any detailed superconduct- 
ing mechanism. While there are problems associated with interpreting various experiments 
within a dx2-y2 context, among the most perplexing from a theoretical viewpoints^ is the 
evident impurity insensitivity of the superconducting transition temperature. Moreover, 
substitution at the rare earth site in YBCO with both magnetic and non-magnetic atoms, 
between the planes of the bi-layer, leads to no variation in Tc, except in the special case of 
Pr. Here we address this fastinating puzzle and demonstrate that for the case of two out of 
phase nodeless s-states in an isolated bi-layer configuration, Anderson's theorem applies to 
inter-plane substitutions. This provides additional motivation for studies of this particular 
pairing state. 

An outline of the paper is as follows. In Section |I| we focus on the bilayer problem 
as a prototype of the iV-layer system. We introduce a generalized model for the pairing 
interaction (Section [11 A] ) and then point out some general properties of solutions which 
arise in these systems (Section [11 B| ). In Section [II Q we examine the consequences of the 
bilayer structure on Josephson tunneling experiments between cuprates and conventional 
superconductors as well as the problem of twinned samples. We then construct (Section [II D[) 
the free energy functional in the Landau-Ginzburg limit, addressing the effects of interlayer 
tunneling and orthorhombicity. Impurity scattering will be discussed in Section [II E| , where 
we investigate the circumstances under which Anderson's theorem holds for the case of two 
out of phase s-states. 



In Section III we turn to the more technical issues of a general weak coupling treatment of 



the problem of two dimensional layers. The reader uninterested in the more mathematical 



details can skip directly to the conclusions in Section |V|. The generalized gap equation in 
the A^ layer problem is derived (Section [III A[ ) and the two competing states are discussed 
in Section [III B[ These correspond to natural analogues of the in-phase and out-of-phase 
states of the bi-layer system. We conclude with a general discussion (Section [III Q) of the 
role of van Hove effects. Two states are found to take maximal advantage of the van Hove 
singularities: these are the d^2_y2 and nodeless s-states. Some of the technical aspects of 
these calculations as well as a complete solution of the bilayer problem are relegated to 
Appendices ^ and respectively. 

Our conclusions are presented in Section |rV|. 

II. THE N = 2 CASE 

A. Gap Equation and Model Interaction 

We initially study the simpler bi-layer system (A^ = 2) in order to establish notation 
and to discuss physical results in a more familiar context. Although only singlet intraband 
pairing is considered in the main part of the text, a more complete weak coupling calculation 
of the bilayer problem will be presented in Appendix p. Many of the interesting features of 
bilayers arise in the general A^-layer problem with only slight modification. However, this 
more general case will be discussed in detail in Section |T|. Section [II B[ focuses primarily 
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on the linearized gap equation. Throughout Section |T|, particular attention will be devoted 
to elucidating the role of orthorhombicity. 

The gap equations for two copper oxides planes are given by 

A, + A.^-I:M,(a.^^ + A-^^) (la) 
A,-^.^-iv.(A,^^^-A.^^^^^ (lb) 



where we define 



El = el + \A^f (2) 



This equation was previously discussed in Ref. Its general derivation is presented in 
the A^-layer context in Section P^II A] . The two order parameter components and A_ 



refer to pairing on the bonding and anti-bonding bands of the Fermi surface respectively. 
Throughout this paper we will transform between the band and layer indices. The two gap 
parameters may be related to order parameters in the layer language by the equations 

A+ = A|| + Ax 

A_ = A|| - Ax ^"^^ 

where Ay refers to singlet pairing of electrons within individual layers and A^ refers to 
singlet pairing of electrons on different layers. Finally, the normal state energy dispersion 
is given by e± = ^ =F t_L, where t± is the interlayer hopping matrix element and ^(q) is the 
single particle dispersion within the copper oxide planes. It should be noted that allowing the 
interlayer tunneling to depend on q in the manner suggested by Chakravarty and Andersonil 
as well as Andersen et al.E2 

^ (cosg^ - cosg^)^ (4) 
has no qualitative effect on the results of this paper. This arises because for the model of Ref. 



28| the inter-plane tunneling deviates most from the constant value t± along the diagonals of 
the Brillouin zone (where Eq. (^) vanishes). It is in precisely this region where the density 
of states contributions are smallest, since, as will be discussed in more detail below, they lie 
away from the four van Hove points. 

As noted above, the bilayer system admits two other pairing states, an interlayer triplet 
state A^ which is odd under interchange of layer labels and a second intralayer singlet state 
Ay where the order parameter has opposite phase on each of the two layers. Both of these 
correspond to inter-band pairing, i.e. pairing between electrons on different sub-bands of the 
Fermi surface. Thus they appear as off-diagonal elements of A in the band representation. 
Because this inter-band pairing is associated with a second phase transition (when t± is 
finite) as well as a triplet state, it appears to be currently of less physical interest than the 
singlet, intra-band pairing which we consider here. For the purposes of completeness, these 
more exotic states are considered in Appendix 

While our analysis of Eq. ([^) is expected to be generally valid, we present here a weak 
coupling treatment, in large part because our emphasis is on the symmetry of the order 
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parameter, rather than on the detailed values of the transition temperatures. Furthermore 
we adopt a model interaction which is believed to be sufficiently generic so as to encompass 
most pairing mechanisms in the literature. Thus it includes the well studied spin fluctuation 
interaction which is stronglvpeaked at quasi-momentum transfers of Q = (vr, 7r),S as well as 
phonon based mechanisms.E^ Alternative scenarios^ may be associated with intermediate 
Q values. In the present model, the in plane and inter-plane pairing interactions are given 
by 

^i(q,q') = A||XQo(q-q') /.n 
i^±(q,q') = AxXQo(q-q') ^ ^ 

where we define the generalized susceptibility, which depends on the momentum transfer 

XQo(q) = T T, TT — } TTT, } TVv^ 

Q=RQo [1 - ^0 (cos(g^ - Q^) + cos{qy - Qy))\ 

This model incorporates several parametrizations, allowing the variation of the peak location 
through Qo, the peak strength via Jo, and the coupling strength and sign through A. We 
consider Jo to be below the critical value 0.5. The interaction contains the underlying 
symmetry of the lattice. This is imposed by summing over all elements R of the point 
group, h being the order of the group. Thus when Qo does not lie on a symmetry element 
of the lattice then the interaction will have eight peaks in a tetragonal ^4/1 lattice and four 
peaks in an orthorhombic lattice. 

Our aim is to incorporate some of the key features of the known band structure in 
various cuprates. Thus the lattice point group will be assumed either to be in the case 
of LSCO or D2h in the case of YBCO and BSCCO. We investigate two different kinds of 
broken tetragonal symmetry in which either principle axes are retained as symmetry planes 
(orthorhombic) in the case of YBCO or the diagonals are retained as planes of symmetry 
(rhombohedral)0 in the case of BSCCO.@ The band structures which we use throughout 



this paper are taken from Ref. In plane energy dispersions of the normal state are given 
by an expansion of the form 

e(q)=Et.r/,(q) (7) 

1=0 

where the basis functions rji are listed in Table |. Also indicated in the table are the various 
band structures we will be considering. It is important to note that the degree to which the 
tetragonal {D^h) symmetry is broken is reflected in the last two parameters of the table which 
we call and ty. The size of these orthorhombic contributions was chosen, for illustrative 
purposes, to be small, but otherwise arbitrary. 



B. Order Parameter Symmetry in Bilayers 

1. Properties of Mono-layer Solutions 

Solutions to the single layer problem provide intuition into those of the A^-layer case, 
particularly when the interlayer hopping is relatively small. It is useful, therefore, to charac- 
terize the allowed order parameters in the one layer limit for each of the three bandstructures 
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discussed above, before moving to the bi-layer case. In the one layer system the functional 
form of the superconducting order parameter depends only on the sign and the shape (as 
distinct from the magnitude) of the pair interaction and on the electronic structure. To 
illustrate these pair interaction and electronic effects, we vary the position of the peak loca- 
tion Q and plot the associated form of the order parameter in this two dimensional space, 
for each of three model bandstructures. The interaction is assumed repulsive (A > 0) in this 
and subsequent phase plots; attractive interactions in general give rise to nodeless s-wave 
solutions for any value of Q. The magnitude of the interaction plays no role in determining 
the symmetry of the order parameter in the single layer problem, it results only in a variation 
of the transition temperature. By contrast, in the two layer case the (relative) magnitudes 
of the in plane and out of plane interactions enter via Eq. (|l|) in an important way. We will 
examine this effect subsequently. 

Other parameters of our model interaction have a lesser influence on the form and regions 
of stability of the various solutions. Variation of the peak width Jq affects Tc to a greater 
extent than it does the symmetry of the actual solution. The component of the interaction 
which most influences the symmetry of the superconducting order parameter is the degree 
to which the vector Q connects positive and negative lobes of the pair wave function on the 
Fermi surface.0 As will be illustrated under more general circumstances below, near half 
filling, when only nearest neighbour hopping is included, the points at which the d^^-y^ state 
has maxima are separated by Q = while the extrema of d^y are separated by (7r,0) 

and the eight lobed Sxy state has positive and negative lobes separated by (vr/2, vr/2). These 
are precisely the points in the space of Q near which these various solutions are favoured. 

Figure |1| represents this "phase diagram" for the LSCO case. It is a map of the irreducible 
representations of the tetragonal, lattice which correspond to the solutions with the 
highest Tc. We stress that the simple functional forms of the order parameter, such as are 
implicit in the notation "c?a.2_y2" etc., can be misleading. The gap parameters may involve 
important contributions from higher order terms in the associated representation. Moreover 
when tetragonal symmetry is broken, states may be highly admixed with states of other 
representations. We will use this notation with care, applying it only to certain solutions of 
systems with tetragonal symmetry. When necessary we will avoid any ambiguity by using 
the common group theoretic nomenclature for the various irreducible representations.B 

Because proximity to the van Hove singularity is found to play such a crucial role, we 
present two limiting cases corresponding to choices for the Fermi energy near (^a) and far 
from (0b) the van Hove singularity. The order parameters shown plot the actual solutions 
for Q = (vr, vr), at which the Big ( "d^^^yi" ) solution is found, at Q = (7r/2, 7r/2), where an 
eight-lobed Aig solution, the so called " s^y" state, exists and at Q = (vr, 0), where the 
"dxy" state is most stable. Moreover, these symmetries persist over an extended region of 
Q-space as outlined by the solid lines in the figures. States with increasing numbers of nodal 
planes, including states in the A2g irreducible representation, appear as Q — (0, 0). Indeed, 
because this region is so complicated (in the case of repulsive interactions), no detailed 
solutions are indicated. The region where the "(ij.2_j^2"-like solution persists is shaded in 
this and subsequent figures. 

Figures (|I]a) and (|l|b) demonstrate that proximity to the van Hove singularity leads to 
a modestly extended region with Big or "dx^-y^" symmetry. Moreover as the Fermi energy 
approaches the singularity, the lobes on the "(i^2_y2" state become sharper. In this way the 
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order parameter takes full advantage of the four van Hove points which coincide with each 
of the four maxima in this gap function. Since the gap function is sharply peaked in the 
direction of the principle axes it is very small in a large region about the nodal lines along 
the diagonals of the Brillouin zone. A large gapless region will have consequences on the 
thermodynamic properties of the samples at low temperatures and may have been observed 
in ARPES measurements on BSCCO.li 

Breaking of the tetragonal symmetry, as in the case of the orthorhombic lattice of YBCO, 
serves to further enhance these van Hove effects. This is shown in Figure ^ where, as in the 
LSCO case above, the gap equation solutions are plotted at the special symmetry points 
and the lines indicate the boundaries of the region over which the representative solution 
persists. Solid lines separate states belonging to different irreducible representations and 
dotted lines separate states of different symmetry within a representation. Figures (^) 
and (||b) demonstrate clearly that the phase space occupied by the "dx^_y2" solution is 
dramatically increased in the vicinity of the van Hove point. Furthermore, the lobes of this 
gap function again become more extended along the principle axes, compared with similar 
solutions farther from the van Hove point. 

Finally, the results in the BSCCO case are presented in Figure ^. We see here that, 
contrary to the two previous cases, the phase space in which the gap has the same sign 
along the a and b directions (i.e, the totally symmetric irreducible representation in the 
lattice) occupies a somewhat larger region of parameter space when the Fermi level is near 
the van Hove point. This is the consequence of a stabilizing admixture of an isotropic s-wave 
component. On the other hand, states with an a, 6-axis vr phase shift dominate the phase 
space away from the van Hove points. 

The above three figures lead to an important conclusion: the "(i^2_j^2" solution exists in a 
large region of parameter space away from Qaf- Therefore this state should be considered as 
an obvious candidate form for the order parameter in any mechanism involving pairing via 
a predominantly repulsive interaction. Stated alternatively, observation of this form for the 
gap in an experiment is in itself not proof of any specific pairing mechanism. Furthermore, 
it also appears from these figures that, quite frequently, a single mechanism can give rise to 
widely varying forms of order parameters depending upon the details of the band structure, 
doping level, lattice symmetry and other subtle features of the material. 



2. Bi-layer Effects 



In bi-layer materials there are two possible states which compete as solutions to the gap 
equation. These may be crudely classified by referring to the states as "inter-plane and 
intra-plane dominated solutions" . These two competing sets of solutions were discussed in a 
previous short communication by our group,! as well elsewhere. By setting t_L = in 

the gap Eq. (|1]) we see that Ay and Aj^ uncouple from each other at Tc. Superconductivity 
can thus be either due to interlayer pairing or due to intralayer pairing depending on the 
relative magnitudes and signs of Ay and A_l. In the band representation in-plane pairing 
yields a solution which has the same phase on both bands while the interlayer pairing 
solution undergoes a vr change in phase from one band to the other. This sign change 
follows directly from Eq. (^) in the extreme limits in which one or the other order parameter 
(in the layer basis) is the larger. 
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It is important to stress that we will refer to these two competing solutions throughout 
this paper. They occur in bi-layer as well as A^-layer systems. In all cases where t± is finite, 
one solution is stable and the other metastable. Thus, although one may refer to their 
associated transition temperatures, only the larger of the two has any physical meaning. 
As a consequence, there is nothing to prevent these two competing states from having a 
different order parameter symmetry. 

As in the one layer case, repulsive interactions are found to lead to nodal states such 
as the (i-states or multi-lobed s-states discussed above, while attractive interactions require 
nodeless states of s symmetry. In the presence of two different types (signs) of in plane 
and out of plane interactions, these effects provide a rich and complex structure for the 
competing order parameter states. Under these conditions and when there is no coherent 
inter-plane hopping, the two solutions may belong to different symmetry groups. When t± 
is finite, the solutions couple while still preserving the underlying lattice symmetry. Thus 
mixing of s- and d-states within the linearized gap requires both orthorhombicity and inter- 
plane hopping. More general s-d mixing effects will occur at higher (quartic) order in a 
Landau- Ginzburg expansion, and will be discussed in more detail in a subsequent section. 

To illustrate the behavior of the bi-layer order parameter and the two competitive states, 
we consider for concreteness the case of an attractive interplane interaction and a repulsive 
intraplane interaction peaked at Q = (vr,7r). This model may be viewed as derived from 
a generalized spin fluctuation model, where the inter-layer or c-axis spin susceptibility is 
included as a perpendicular component of the Dairing interaction. This case was previously 
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discussed by a number of different groups.&clOO'ES While we choose this peak location simply 
for illustrative reasons, our results can easily be extended to other values of Q. Under these 
circumstances, the order parameter in the band picture has the schematic form 

(8) 

If the lattice is tetragonal then either s± or d\\ will vanish. Thus when A_l is zero the solution 
will have (i-symmetry. Moreover, both band gaps are in phase. Similarly when Ay vanishes 
the symmetry is s-like and both gaps are out of phase. The transition from one type of 
solution to the other is governed by size of X±/X\\- Frequently the cross over from the in- 
plane to the inter-plane dominated regime may occur despite the fact that the ratio A_l/A|| 
is small compared to one. What determines this cross-over is the nature and sign of the 
competing interactions and the bandstructure. 

If tetragonal symmetry is broken even more complicated situations obtain. Both com- 
ponents can be nonzero simultaneously and the transition is a smooth one as a function 
of increasing A_l/A||. The regions where the various types of solution exist in YBCO are 
illustrated in Figure | for the example discussed above. Plotted on the vertical axis is the 
fractional hole doping x. We estimate the physical values of this parameter, for a range of 
different oxygen stoichiometrics in YBCO to vary from 0.1 to 0.3. The flgure illustrates the 
parameter regimes where the dominant component of the order parameter on the two bands 
is {d, d), (s, d) and (s, — The region below the solid line involves states in which the sum 
of the pair wave functions on the two bands changes sign under 7r/2 rotation of the axes and 
can thus yield vr phase shifts in a, 6-axis Josephson tunneling experiments.EirO Moreover, 
the shaded regions denote states in which the sum of the two order parameters changes sign 
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even though one or both are nodeless. This issue will be addressed in detail in Section [II C 



The inset plots the density of states as a function of energy. There are actually four van 
Hove points associated with orthorhombic splitting of the two bands. For clarity the doping 
level at which the first of these occurs is plotted in the main portion of the figure as a dotted 
line. Note that the region labelled by (s, — s) occupies the largest fraction of phase space 
when the band is near the van Hove points in general and in particular when one sub-band 
of the Fermi level is closed about the F point of the reciprocal lattice and the other about 
X. 

The evolution of these solutions with increasing A_l/ \\\ is illustrated in Figure ^. Indicated 
here is the calculated shape of the gap functions for the two bands of the bi-layer. Their 
relative phase is also noted. It can be seen that the pair wave function on each band has a 
greater amplitude along one symmetry axis than along the other, and the axis along which 
the pair wave function has its maximum differs from one sub-band to the other. These 
results are plotted for very small orthorhombicity and for Fermi energies close to the van 
Hove point. On the basis of the latter assumption, the small orthorhombic effects are greatly 
magnified. The figure labels (a)-(d) correspond to (a) two in phase c?-wave dominated states, 
(b) one band nodeless and the other nodal, (c) two out of phase nodeless states which show 
greatly elongated lobes, and (d) to two out of phase s-wave dominated states which are 
more isotropic. The relative phase space which these solutions occupy can be seen from the 
previous figure. 

Moreover, the detailed shape of the order parameter has important physical conse- 
quences. Thus for case (c), for example, the solutions are in the inter-plane regime where 
A± is dominant. Here, however, the order parameter is nodeless and yet the combination 
of the two gaps has features of a (i-wave solution in a, 6-axis as well as c-axis Josephson 
tunneling experiments.^ Case (b) is associated with one s-like and one d-like state. It will 
exhibit both power law behavior in thermodynamical properties, as well as the vr phase shifts 
in corner junction experiments, which are also seen in cases (a) and (c). These Josephson 
experiments and their relation to the order parameter symmetry will be discussed in more 
detail in the next section. 



C. Josephson Tunneling: Phase Coherence Across Domain Boundaries 

Josephson tunneling experiments have been key to elucidating the order parameter sym- 
metry in the YBCO cuprate family.&i3@ Within this class, two types of measurements 
have been performed. These involve a, h- and c-axis tunnel junctions. In the former cate- 
gory SQUID geometries have been investigated.lii^ These consist of two junctions with Pb 
counterelectrodes, whose interference pattern yields information about the relative phase of 
the order parameter along various axes in the a, 6-plane. In "corner junctions", vr phase 
shifts of the order parameter between the a and h axes have been observed and attributed 
to the sign change upon a 7r/2 rotation of the wave vector, within a dx2-y2 state. This 
interpretation is further confirmed by "edge" junctions which probe the interference along 
a single face of the material; here no phase shifts are found to be present. These conclu- 
sions are the same for both twinned and untwinned crystals. A variant on this geometry 
are the ring experiments^ consisting of YBCO segments with different grain boundary ori- 
entations. Observation of a non-zero, half integer spontaneous flux threading the ring, for 
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specific orientations of tlie grain boundaries, again provides support for a d^2_y2 state. 

An alternate Josephson tunneling geometry has been investigated by Sun et al@ These 
authors observe finite c-axis Josephson currents between a Pb counterelectrode and a YBCO 
sample. For the simplest representation of the d^i-yZ symrnetry, no Josephson current should 
be present. While it is possible to invoke orthorhombicity^a to explain these non- vanishing 
results, this hypothesis seems unlikely since the magnitude of the current is evidently not 
sensitive to whether twins are present. 

To the extent that twinned YBCO can be treated as a tetragonal system, these two types 
of experiments appear manifestly incompatible, tt phase shifts are associated with a sign 
change of the order parameter in the a, b- plane. The c-axis current, which reflects an average 
of the order parameter in the a, 6-plane, must necessarily vanish. These observations are 
not altered when bi-layer effects are incorporated. While there has been some attentionMi^l 
paid to the role of chains in addressing the data, it should be noted that recent experiments 
on reduced oxygen (where the chains are expected to be highly interrupted and, 

therefore, irrelevant to the transport and superconductivity), seem to reproduce the same 
behavior as in the optimally doped material. In summary, if both a, b- and c-axis Josephson 
measurements prove to be correct, any resolution of this issue will probably revolve around 
a deeper understanding of twinning effects. 

In a Josephson experiment, the measured Josephson current is the sum of the currents 
established between each band of the cuprate and the superconductor to which it is coupled 

J=J+ + J- (9) 

where each component has the usual Ambegaokar-Baratoff form0 

= AiA^^E N + + A^^)]"^ (10) 

The resistances, r±, can be calculated by determining the single particle tunneling matrix 
elements^ 

r,%, = JdT{<Pi^\j{r)\<p§) (11) 

where the integration is carried out over the intermediate region. The wave-functions 0^ and 
(j)^ are the single particle wave-functions from the left and right hand side of the junction 
respectively, and j(r) is the usual current operator. The fact that the electrons of the 
cuprates are localized on the copper-oxide planes requires that they be described by wave 
packets with some finite spread in momentum along the c- axis direction. A microscopic 
treatment of the tunneling process is complicated by a number of issues, which we cannot 
address here: the propagation of electrons out of a region of localization in the c-direction 
will lead to scattering effects as the tunneling pair enters a more isotropic material. These 
may influence the tunneling from both bands to a substantial degree. (It should also be 
noted that the order parameter at the surface may be modified from its bulk form.ii) For 
definiteness, we first consider the case of coupling to a conventional superconductor such as 
Pb, and ignore these complications. 

If we make the simplifying assumption that the Josephson coupling is the same for 
both the symmetric (+) and antisymmetric (— ) bands, it follows that tt phase shifts in a 
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SQUID experiment will be observed in a substantial region of the phase space, corresponding 
to all bi-layer parameter sets which lie below the solid line of Figure Moreover, for a 
significant fraction of these (particularly in the vicinity of the van Hove points), the states 
in question differ from two in phase ci-like states. They are described by (s, d) or (s, —s) 
dominant combinations. In this way, the measurement of tt phase shifts cannot be uniquely 
associated with d states in a multi-layer system, provided orthorhombicity is also present. 
This observation makes it all the more important to repeat the various SQUID experiments 
for either mono-layer or tetragonal materials. 

In reality, there is some asymmetry in the Josephson coupling to the symmetric and 
anti-symmetric bands. This asymmetry derives from (1) density of states effects related to 
the van Hove points and (2) wave- function structure: wave functions associated with the 
anti-bonding band of the bilayer cuprate have opposite phase on the two layers and therefore 
a nodal plane exists between the layers. The first effect tends to enhance Josephson coupling 
preferentially to the antibonding and the second to the bonding band. The net contribution 
cannot be calculated with any certainty, but it is reasonable to conclude that the Josephson 
coupling is appreciable for both bands of the bi-layer system. 

These remarks can be addressed in somewhat more detail. In most of the high cuprates 
for which band structure data is available, the Fermi surfaces are closed about the (vr, vr) 
point. As a consequence the anti-bonding band lies closer to the van Hove points and thus 
has a higher density of states. This density of states contribution will then lead to a relatively 
stronger Josephson coupling for the anti-bonding or antisymmetric band. On the other hand, 
the amplitude of the pair wave function is expected to be greater on the even symmetry 
bandii and this effect may lead to a greater contribution to the Josephson current than 
for the even band. Moreover, it has been claimedB on the basis of symmetry arguments, 
that tunneling from the anti-bonding band into an s-wave superconductor will lead to an 
appreciable reduction in the current contribution from this band. However, some Josephson 
coupling is expected to remain; the matrix elements TfTj./ only vanish under the special 
circumstances, when the centers of the two wavefunctions of Eq. (|ll]) co-incide throughout 
the boundary region. In general, this matrix element will be nonzero, if for no other reason 
than because the lattice constants of the two materials are unlikely to be commensurate. 
The situation is depicted in Figure ||. 

Josephson tunneling from one YBCO crystal to another, as well as across twin bound- 
aries, is even less amenable to microscopic theory. In many respects twinned materials be- 
have like single crystals, with similar transition temperatures, thermodynamics and Joseph- 
son currents in corner and c-axis junctions. It is clear that little is understood at a detailed 
theoretical level about the nature of the twin boundary. What seems to be less ambiguous, 
however, is that there is little if anv pinning of the critical current (in low magnetic fields) as 
it flows between twin boundaries.E3 This would suggest that, whatever the order parameter 
symmetry, the phases tend to line up with -|- and — lobes adjacent. 

In the context of a bi-layer system, this situation is more complex, since the phases 
are associated with multiple bands. In Figure |^ we schematically plot the two alternative 
scenarios for the (s, —s) (orthorhombic) bands (parts (b) and (c)), as well as the generally 
expected behavior for the {d, d) ). Of the two scenarios (b) and (c), only the latter 

would preserve the vr phase shift behavior across a twin boundary. This scenario would also 
be compatible with a low twin boundary pinning of the critical current. 
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Despite these phenomenological arguments, case (b) cannot be ruled out on microscopic 
grounds. In the extreme hmit where the boundary can be treated as a conventional SIS 
junction, the locking of the phase of the order parameter across the junction occurs via 
Josephson coupling. If the bilayers on either side of the junction are properly aligned it 
might be expected that only bands of equal symmetry couple together since the matrix 
element in Eq. (|ll]) vanishes otherwise. Thus the even bands on either side of the junction 
would couple together as would the odd bands. As a consequence the positive lobes of a 
(A+, A„) = {d,d) solution would point in the same spatial direction on both sides of the 
twin boundary (scenario (a)). On the other hand the (s, —s) state would have the component 
of positive phase in one twin aligned with the component with negative phase in the other 
twin (scenario (b)). This would lead to a substantially reduced Josephson current at a 
macroscopic junction due to averaging over the twins. If instead of an SIS model, one 
argued that the order parameter lobes were required to vary in the most continuous fashion, 
one might conclude that the dominant component should also be continuous. This, too, 
would lead to scenario (b) for the (s, —s) states and (a) for the {d, d) configuration. 

Nevertheless, it is also possible to assume that interlayer pair breaking effects become 
considerable in the boundary region and that Ay retains its coherence between twins so that 
scenario (c) obtains. Moreover, stacking faults, lattice defects and other complexities can 
invalidate any of the above simple models of twinned crystals. In summary, the nature of the 
order parameter variation across a twin boundary is quite complicated in one layer systems 
and sufficiently complex in the bi-layer case so that no clear conclusions can be drawn at 
this time. 



Thus far we have investigated only the linearized gap equations, which are necessarily 
restricted to the vicinity of the transition temperature. Additional effects may occur below 
Tc associated with the transition to states with other order parameter symmetries. While 
there does not appear to be experimental evidence for additional phase transitions, there is 
considerable information contained in studying the more general situation. In this section 
we derive the appropriate Landau-Ginzburg free energy functional in terms of both the layer 
and band indices. The behavior at quadratic order gives further insight into, and serves to 
validate the results discussed in the previous sections. There have been several discussions 
in the literatureBS of phenomenological forms for the bi-layer free energy. Here we proceed 
from a microscopic basis. Note that while this discussion refers to a bilayer structure, it is 
readily extended to the case of general N following the results of Section [111 Ij| . 

In the case under consideration (when only pairing of electrons on individual sub-bands 
of the Fermi surface is considered), the superconducting state of the bi-layer is described by 
a two component order parameter : A_|_ and A_. The free energy, Fg, of the superconducting 
state then consists of the sum of terms 



where Fs^± is the contribution to the free energy of an individual sub-band. In the Landau- 
Ginzburg limit the free energy difference between the normal and superconducting states is 
given by 



D. Landau-Ginzburg Free Energy Functional 



Fs 



Fs,+ + Fs, 



(12) 
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T /• dS , , ,2 0.0533 f dS , , ,4 

Tc Je±=Ep [ZTTj-^VF [KbTcY Je±=Ep [ZTTj-^VF 

It should be stressed that, although the two sub-bands contribute independently to the free 
energy, the two order parameters are coupled via the gap equation. Both become finite at a 
common Tc. 

By use of the gap equations the the Landau-Ginzburg free energy functional can be 
recast into a form which contains a quadratic part 



a+|A+|' + a„ |A_r + S [A+A*_ + c.c.j (14) 
and a quartic contribution 

/3+|A+|^ + /3_|A_|^ + (A+A*. + ex.) f/i+|A+|2 + /i_|A_p) (15) 



Here we have assumed that each component of the order parameter can be separated into a 
complex magnitude times a normalized real function over the Fermi surface: 

A±(q) = A±^±(q) (16) 

The magnitude A±. appears explicitly in the Landau-Ginzburg free energy expansion while 
'ip± determines the co-efficients in the expansion. The co-efficients 6 and fi±. vanish if A_|_ 
and A_ belong to different irreducible representations of the lattice point group. 

Thesecoefficients have been discussed by Varma in the context of a particular pairing 
scenario.E2l While our free energy contains the same class of terms as that presented in Ref. HO, 



because we have explicitly removed the wave vector dependences via Eq. (p!6[) above, it is 
not straightforward to determine the conditions for order parameter sign changes under a 
tt/2 rotation of the lattice. Nevertheless, the relative phase of the two order parameters 
can be determined variationally, by minimizing the free energy. If the nontrivial coupling 
parameters {6 and fi±) are finite then the phase difference between A+ and A_ can be at 
most or tc. This last result is consistent with the discussion of Section 111 B 21 and the 



appropriate sign depends on whether in-plane or interplane correlations are dominant .S 

It is useful to transform the above free energy functional to the layer representation 
using Eq. (|^). The functional is expressed in terms of Ay and A^ and two explicit transition 
temperatures can then be associated with the intra- and inter-plane interactions. In this 
context the mixing between various symmetries can be understood in a more direct way. In 
this representation the free energy is given by 



a|||A|||2 + ai_ \A^f + 6' (A||A1 + c.c.) + Ay 1^ + /^^lA^I^ + /3'|A|||2 \A^f 

/(AjjAf + C.C.) + (A||A1 + c.c.) I A|,|2 + ;xx I A^''^^ ^ ^ 



-7 



Here the coupling terms 6' and vanish if Ay and Aj^ belong to different irreducible 
representations of the lattice symmetry group. Note that one important effect of orthorhom- 
bicity is to require that the coupling parameters 6' and yU|j/± change sign upon the interchange 
of the a and b axes of the crystal. These same terms also vanish in the limit as t± goes to 
zero. 

By direct calculation from the microscopic theory we find that the quartic cross-term 
7' is positively so that the phase difference between Am and A_|_ is zero (if 6' < 0) or vr (if 
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5' > 0).E3 Studies of a related Landau-Ginzburg free energy functional have been presented 
by Kobuki and who discussed the mixing between a d-wave Ay and an s-wave A_|_ in 

an RVB based theory. In their approach, mixing between these two components was brought 
about by a self-consistently determined orthorhombic strain. This led to the introduction of 
additional terms into the free energy which could result in a negative 7'. As a consequence 
an s + id state was produced. In the present theory, in contrast, mixing between A_|_ and Ay 
occurs at quadratic rather than quartic order, as a consequence of finite tj_. Furthermore, 
in the absence of tetragonal symmetry breaking only s + d mixing occurs (at quartic order). 
Additional arguments against an s + id state were presented by Normand et al.@ 



E. Impurity Effects 

In this section we discuss the nature of impurity or pair-breaking effects in bi- layer 
systems. At the heart of this issue is the paradoxical observation that all substitutions at 
the rare earth site, which sits between the bi-layers (except for Pr), leave unaffected. 
Rare earth substitutions with or without local moments and in a disordered or ordered form 
make no difference to the superconducting transition temperature. Previously it has been 
argued that the even more general insensitivity of the cuprates to impurity substitution 
is incompatible with anisotropic or d-wave superconductivity.!!^ Here we investigate the 
complexity introduced into this problem by the presence of a bi- layer order parameter. 
Since it is relatively straightforward to generalize to the magnetic case, for definiteness, we 
concentrate on the case of non- magnetic impurities. 

A new aspect of the present work is the consideration of an isolated bi-layer, rather than 
a coherent stack of bi-layers. Our starting point is a necessary first step in a treatment 
of incoherent coupling along the c-axis. As a result of this assumption the configuration 
averaging process (which restores the underlying translational symmetry of the lattice) is 
different from that discussed in Refs. |T3i^ As has been noted elsewhere,^ there are two 
types of impurities which must be considered: intra and inter-layer substitutions. For an 
isolated bi-layer, inter-layer impurities represent the more interesting case, since the relevant 
equations decouple when written in the band basis. Thus a state such as an isotropic 
(s, —s) state will be insensitive to inter-layer impurities (as will the in-phase (s, s) state). 
By contrast, intra-layer impurity effects involve processes which couple the two bands. In 
this way, the (s, —s) state exhibits intra-layer pair-breaking.lli Moreover, the same concerns 
that were raised earlierSS about the d-wave order parameter apply to the bi-layer case with 
either type of impurity. Thus the (s, — s) state emerges as the leading (non-trivial) candidate 
state for resolving the paradox concerning rare earth (non-magnetic) substitutions in the 
cuprates. 

We begin with the standard treatment of scatterers within individual planes using the 
Born approximation.^ The impurity Hamiltonian has the form 

^Sane = EE«Kq)4,.(P + q)c..,.(p)e-^''-^- (18) 
m q,(T 

where m labels impurities located at Rm on layer im and cr is a spin index. The impurity 
self-energy is given by 
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S^,r,r'(p) = -^"-imp I"^(P ~ Q) 1^ ^a;,fc" (p)^r,r' (19) 
r" ,q 

where Q is the averaged Green's function and riimp the concentration of impurities per layer. 
We assume that the gap in the absence of impurities is given by A-i-(q) = A±'?/'±(q) where 
are functions normahzed appropriately over the corresponding sub-band of the Fermi 
surface and that the impurity renormalized gap is given by A-i-tlj±{q). The self consistent 
equations then become 

1 LJl 

^i = ^i + TZ2^ I. ~, (20a) 



r' 



where r is the usual scattering time and 

Here is a combined measure of the anistropy of the order parameteiH and impurity po- 
tential. This coupling constant varies from to 1. The latter is appropriate to the case of 
a totally isotropic order parameter. 

It follows from Eqs. (EUp that intra-layer impurities enter the coupled self consistent 



equations via a mixture of the two band contributions. In this way they lead to pairbreaking 
in all instances, except for the special case of two in phase, isotropic s-states. Thus the (s, —s) 
states experience a reduced Tc in the presence of these impurities.0 

We next consider inter-plane scatterers. For the sake of generality these impurities are 
assumed to scatter electrons within as well as between planes. The impurity Hamiltonian is 
given by 

^Ltoplane = E E [^11 (^Uci,. + 4,.C2,.) + (cl^C^,. + C.C.)] e-'^^-^'" (22) 
m q,(T 

It is important to note that this Hamiltonian is diagonal in the band language and that no 
average is taken over sites in the vertical direction. Thus the renormalized self-energy and 
order parameter satisfy equations which are decoupled in the band index 

0Ji,± = uji + , (23a) 

A± = A± + , = (23b) 

2r,y-2;^ 

Here care must be taken to preserve the band labels on r and g. This decoupling of the two 
bands leads to the conclusion, stated earlier, that isotropic order parameter sets (whether 
in or out of phase) experience no pairbreaking from inter-layer impurities.!! 
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In summary, we see that one way of avoiding the strong pair-breaking generally associated 
with substitutions at the rare earth site between the bi-layers is to consider states with the 
symmetry (s, — s). For this reason, along with the tt junction behavior discussed in the 
previous section, this state should be considered as a potentially interesting candidate for 
the YBCO system. 



III. FORMULATION AND ANALYSIS OF THE iV-LAYER PROBLEM 

A. Gap Equation in the Band Representation 

In this section, we treat the general A^-layer problem. The system under consideration 
corresponds to a stack of decoupled A^-layer structural units, each layer of which consists 
of a two dimensional copper oxide plane. In the limit of infinite N, we recover essentially 
the usual Bloch wave description of a collection of copper oxide mono-layers, aligned along 
the c-axis.0lli Just as in the bi-layer case, we find that there are two competing states. 
Depending on the relative size of the intra- and inter-plane coupling constants; one of the 
two is stable, while the other is metastable. These two limits correspond to inter- and intra- 
layer dominated states. In the band language the latter are in phase and the former out of 
phase. Because the system contains such bands, the inter-layer gap parameters should 
be viewed as sinusoidally modulated with varying band index, as will be illustrated in more 
detail below. 

Just as in the bi-layer case, the charge carriers within individual planes and on adjacent 
planes within a unit cell are assumed to interact via a non- retarded pair potential. Hopping 
of quasiparticles between adjacent planes within a unit cell is determined by the hopping 
matrix element t±. No hopping is allowed between unit cells, as a consequence of our 
assumption that the c-axis coupling is incoherent. As before, only singlet intraband pair 
states are considered. Many of the detailed derivations in this section may be found in 
Appendix |A 1. 



The non-interacting Hamiltonian has the form 



N N-1 



+ C.C.) 



o" Li=l i=l 



(24) 



where i is a layer and cr a spin index. The superconducting order parameter is defined by 

A..(q) = E^^.(q.q')i^..(q') (25) 

q' 

where Vij is the interaction between electrons on layer i and layer j and the anomalous 
Green's function Fij{q) = (cj^(q)cjj(— q)) is antisymmetrized with respect to spin indices. 

The anomalous component of the superconducting Green's function is obtained by per- 
forming the usual matrix inversion 



LUJi — Hn 

-A* LUJi 




(26) 
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The entries in this matrix problem are 2N x 2N matrices with two spin degrees of freedom 
and N layer indices. 

We diagonalize Eq. (^) in the spin degrees of freedom by multiplying the entire equation 
from the right by a2- The equations are then transformed to the band picture by diagonal- 
izing the normal state Hamiltonian ifjv- The resulting energy dispersion in band r of the 
Hamiltonian is given by 

er(q) = e(q) + 2tx cos(^^^j r = 1, . . . , iV (27) 

We will consider only pairing of electrons within individual bands so that the order parameter 
is diagonal in the band representation with components A^. With a fully diagonalized 
Green's function the gap equation for the iV-layer system can now be readily obtained. 
Upon performing the sum over Matsubara frequencies, we find 



A.(q) = -EE^^'^' (q,q')^^tanh 

q' r'=l 2hri{q.) 



(28) 

The energy dispersion of the elementary excitations is given for each band r by the usual 



relation Er = ye^ + |Aj.| . Equation (|28|) is the central equation of this section. In Ap- 
pendix |A3| we generalize this result further by extending it to an infinite stack of layers, 
corresponding to a fully three dimensional lattice. 



B. Solutions of the Gap Equation for Small t± 

In this sub-section we establish the nature of the two competing states which are the 
stable and metastable solutions to the gap equations derived from Eq. (pS]). These are most 
readily introduced by considering first the limit of small t±. For arbitrarily small t±, the two 
solutions become independent and appear with different onset or transition temperatures. 
These two states are respectively associated with pure intra- and pure interlayer pairing. 
Moreover, in this limit analytical results can be obtained, while the more general case of 
non-zero t± is treated numerically. 



After some algebra, which is outlined in Appendix |A 2| , it follows that the solution to the 
gap equation is given by 

A,(q) = A||,o ^ii(q) - 2cos(^^^^j A^.o ^±(q) (29) 

where the two Fermi surface functions ^/'||(q) and V'±(q) satisfy 

r dS 

^11 ^11 (q) = - y Ml q') ^11 (q') (30a) 

r dS 

n^M^) = - /tt^-^ V±(q,q')^±(q') (30b) 

Here the integrations are over the degenerate bands of the Fermi surface and the related c- 
numbers, fl\\ and fl_i, are related to the respective transition temperatures T^\/±fi determined 
by the equation 
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Jl^] . J- (31) 

Finally, the two competing states are associated with taking either one of the two parameters 
in A|| or A^^o in Eq- ( P^ to be zero. 

The above results can be generalized to the case of finite t±, using numerical techniques. 
When t±_ is finite, the two parameters Ay^o and A_l.o can be simultaneously non-zero in which 
case both and %l)± will belong to the same irreducible representation of the lattice point 
group. 

To illustrate these results, we plot the amphtudes on the various bands of the competing 
inter- and intralayer pair states for = 4 and = 7 in Figure ^ As in the bi-layer case, 
if the dominant interaction is attractive then the order parameter will be nodeless, while a 
repulsive interaction will yield a nodal solution with exact form determined by the details 
of the interaction and the band structure (as discussed in Section [TlB 1| ). In panels (a) and 



(b) of Figure ^ the magnitude of the order parameter is plotted as a function of band index 
r for the intraband (a) and interband (b) states for a four layer system. The dotted line 
indicates the analytical solution for = and the histogram bars illustrate the numerical 
results for moderate (comparable in magniude to the separation of the Fermi level from 
the van Hove points). In this way some deviation from the analytically obtained curve is 
seen as the solid bars differ slightly from the dotted line. It is clear that the two competing 
solutions represent a natural generalization of the bi-layer results to an A^-layer system. 
Similar results are plotted for the seven layer system in panels (c) and (d). It follows from 
the figures and the above discussion that in these higher A^ systems, even more complex 
behavior can be obtained, with a range of signs and magnitudes of the order parameters 
associated with the different bands. 

To make this complexity even more explicit, we have considered the case of A^ = 3 for 
the case of a dominant attractive interlayer interaction and a weak in-plane repulsion both 
peaked at (vr,7r). In this case the hopping is slightly larger than in the previous figure. 
Figure | shows our solution to the gap equation and the normal state bandstructure (inset) 
for this 3 band model of orthorhombic YBCO with the two lower energy bands closed about 
X and the highest energy band open. This last band can be viewed as simulating the chain 
band in YBCO,il since conduction in this band is only possible along one principle axis. 
Solution of the 3 layer gap equation clearly shows the mixing between the two components 
of the solution. The middle (plane-like) band has a pure d-wave solution, since interlayer 
pairing contributes very little to the gap on this band, while the other two bands have s-wave 
symmetry. The physics of this 3 band model is equally complex. Because of the dominance 
of an s-wave order parameter component, one expects that the magnitude of T^. is only 
mildly affected by impurity scatterers. On the other hand, a nodal solution on one of the 
bands will yield power law dependences in thermodynamic functions at low temperatures. 
While the above model should not be viewed as a detailed representation of YBCO, it serves 
to illustrate the rich array of phenomena which are associated with multiband systems. 
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C. Basis Functions and Van Hove Effects 



Band structure effects have played an important role in our analysis, particularly when 
the Fermi energy lies in the vicinity of the van Hove singularities. We have seen that 
these singularities distort the shape of the order parameter. They also play a key role in 
determining the relative stability of various solutions to the gap equations. In this sub- 
section we show that, of all the different gap symmetries, two are able to take maximal 
advantage of the van Hove points: these are the nodeless s if the interaction is attractive 
and the "dx2_y2" statesH in the case of a repulsive interaction. It should be noted that there 
is a considerable literature on the effect of the van Hove singularities on raising In this 
paper we focus on the interplay of the van Hove singularity and order parameter symmetry. 
Earlier work^S has shown that in the more general strong coupling picture, states far from 
the Fermi energy may wash out, to some degree, the effectiveness of the van Hove singularity 
in raising T^. While here we use a weak coupling approach to address the order parameter 
symmetry the same qualitative behaviour can be expected to follow in more general strong 
coupling calculations. 

To quantify the van Hove effects we study the linearized form of Equation (^8|) 



where the integrations are performed over segments of the Fermi surface corresponding to 
the different bands and the eigenvalue fl defines the BCS transition temperature 



The Fermi velocity on band r is Vp, and Ep is the Fermi energy. 

We define a complete set of orthonormal basis functions over the Fermi surface,^ which 
are non-zero over only a single band r and assume further that these 'ipi^^\c{) belong to an 
irreducible representation F of the lattice point group. These basis functions satisfy 



The pair wave function and interaction potentialEll are expanded in terms of these Fermi 
surface harmonics as 




(32) 




(33) 




(34) 



(35) 



and 




(36) 



The gap equation is thus reduced to the simple set of eigenvalue problems 



(37) 
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By working in the space of functions defined by Eq. (Q) it is clear that the basis functions 
ipl^'^^ are weighted by the inverse of ^/vp and so regions along the directions of the van Hove 
points give a correspondingly greater contribution in the gap equation. If the interaction 
V^'"^ (q, q') coupling two points on the Fermi surface is repulsive (attractive) then states with 
opposite (same) phase at these two points will be favoured. Thus we can conclude that the 
dx2_y2 basis function will benefit most from the van Hove points for repulsive interactions 
while the nodeless s function will be most enhanced by an attractive interaction. 

The effect of the van Hove singularities in the single particle density of states on the 
superconducting order parameter can be quantified through the pairing density of states 
(PDOS)i 

|A(q)r 



Pa(E,) = (38) 

dS|A(q)r 

The integrals are taken over all bands of the Fermi surface. To illustrate this function and 
its relation to the van Hove singularities, in Figure (p!0|a) we plot Pa as a function of the 
Fermi energy for various solutions to the gap equation. Here we focus on the one layer case 
for clarity. It can be seen that a peak appears at the van Hove point for all irreducible 
representations. It is, however, more appropriate to normalize by the single particle density 
of states Pa=i- This gives a more accurate indication of the degree to which the pair wave 
function is stabilized by the density of states. The result for the four different irreducible 
representations of the D^h lattice of LSCO is plotted against Fermi energy in Figure ([TO|b). 
Here we select that Q which results in a maximal Tc for each representation. In this way, 
we see that the Big solution, ie. the "d^s^y^" state, and the nodeless Aig "s"-wave solution 
indeed benefit much more from the van Hove singularity than do all other states. 

In summary, there are two states which take maximal advantage of the van Hove points. 
These are the nodeless s-wave state which occurs only for attractive interactions and the 
dx2-y2 state, appropriate to the case of repulsive interactions. It should be stressed that in 
general the order parameter will not have these simple functional forms corresponding to a 
single basis function. This is is all the more striking as the Fermi energy approaches the van 
Hove singularity where admixtures of higher order basis functions are most evident. While 



in the presence of multiple bands, the results plotted in Figure |T0| become more complicated, 
the essential features still remain. 



IV. CONCLUSION 

The most important issues in the field of high T^. superconductivity involve determina- 
tions of the order parameter symmetry and the superconducting pairing mechanism. While 
there are, clearly, no definitive answers to be had at this time, this paper has been directed 
towards addressing these two issues. We have emphacized the role of multi-layer effects in 
the cuprates, in large part because the most well characterized material, YBCO, has two 
copper oxide planes. This complexity leads to complications in inferring the order param- 
eter symmetry from various experimental tests. It also suggests that there are different 
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(inter- and intra-plane) channels which should be considered in any microscopic theory of 
the pairing. 

In reference to the order parameter symmetry, we have found that a multi-layer system 
should be characterized by distinct gaps appropriate to each of the multiple bands. A bi- 
layer material such as YBCO has two gaps, a tri-layer, three, etc. In the presence of even 
a very small amount of orthorhombicity, one of the gaps can be predominantly of d, while 
another of s symmetry. (Throughout this paper we refer to d and s states as those which are 
odd or even, respectively, under a 7r/2 rotation of the wave- vector) . One may be node-less 
while the other has nodes. The multiple gaps can be in or out of phase. In this way the 
observation of power law behavior in thermodynamics may reflect on only one of the order 
parameters in question. The observation of vr phase shifts in Joscphson corner junction 
experiments on YBCO must be viewed more widely in this multi- band context. Indeed, we 
have found that this behavior can be associated with two out of phase s-states in the presence 
of (weak) orthorhombicity. This orthorhombicity leads to a strong asymmetry of the s and 
—s states, so that one gap function is elongated along the a and the other along the b axes of 
the crystal. The net Josephson current behaves rather similarly to a d^^-y^ state, although 
the thermodynamical behavior need not exhibit the power laws of this state. Finally, the 
behavior of multi-layer systems in the presence of impurities is similarly complex. Intra- 
and inter-layer impurities suppress Tc in a different fashion. The s, —s state is of interest 
because it obeys an Anderson theorem with respect to inter-layer substitutions. This may 
help explain why substitutions at the rare earth site in YBCO make little or no difference 
to the magnitude of Tc. 

In the process of investigating very generic model interactions for the superconductiv- 
ity, we have inferred information about microscopic constraints on the pairing mechanism. 
dx2_y2-like states are found to be general solutions to the gap equation for repulsive inter- 
actions, in large part because they possess the fewest number of nodes and thereby the 
highest transition temperatures. In this way, they should not be specifically associated with 
a spin fluctuation driven pairing mechanism. Moreover, van Hove effects act to stabilize 
some order parameters over others. Of these the dx2-y2-like symmetry is, again, the most 
notable. Orthorhombicity further enhances this stabilization. Thus for a variety of reasons, 
this state emerges as a natural solution to the gap equation (s) in the presence of repulsive 
interactions. 

While we have emphacized the bi-layer {N = 2) case, we also presented general multi- 
layer calculations which view the c-axis as consisting of decoupled structural units, each of 
which contains copper oxide layers. By contrast, within the unit cell the intra- and inter- 
plane hopping is appreciable and plays an important role in giving rise to A^ distinct bands. 
These N > 2 calculations may be particularly relevant in the context of the Hg and Bi based 
cuprates. In treating the superconductivity, we have included intra- and inter-plane pairing 
interactions in parallel with the above intra- and inter-plane hopping. We demonstrated that, 
regardless of the number of layers A" in the unit cell, there are always two competing states: 
one of which is intra-plane dominated, so that the resulting A^ band- gaps are in phase, 
and one of which is inter-plane dominated, so that the A^ gaps are sinusoidally modulated. 
Small changes in the parameterizations can lead to a transition from one of these states to 
another. Thus it may be inferred that the order parameter symmetry is potentially variable 
from one cuprate to another and from one stoichiometry to another. 
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While the inclusion of these multilayer effects has been seen to introduce considerable 
complexity into the classification of the order parameter symmetry, this complexity is in- 
escapable. As long as the layers communicate via one or two body processes (i.e., via 
hopping or pairing interactions), superconductivity in the high Tc cuprates must include 
these multilayer effects. 

Note Added: After this manuscript was completed we learned of recent experiments 
from C. Tsuei and co-workers in which vr phase shifts have been reported for a one layer, 
tetragonal Tl compound. 
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APPENDIX A: CALCULATION OF GAP EQUATION FOR iV-LAYERS 



1. Derivation 



To begin we need to calculate the superconducting Green's function defined by the 
2N X 2N matrix 



,-1 ^ ( ^i- Hn -A 
' ~ I -A* LuJi + Hn 



(Al) 



The components and A are defined in Section [ill A] . First we diagonalize the normal 
state Hamiltonian. For this purpose we define the set of characteristic polynomials 

DN{i-e) = dei{HM-e^N) (A2) 
and observe that they satisfy the recursion relation 

Di{x) = xDi_i{x) - t]_Di_2{x) 

Do{x) = (A3) 
Di{x) = X 

By using the known properties of the Chebyshev polynomials it is straightforward to show 
that these characteristic polynomials have the form 

Dj,{x) = i-t^f + 
^ ^ ^ ^ sm(g) (A4) 

X = —2t±_ cos(g) 

From the zeroes of these polynomials we can determine the normal state energies 

eM) = e(q) + 2tx cos(^^^ j r = 1, . . . , iV (A5) 
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An orthonormal set of eigenstates of H^^ can be found similarly. The i*'* component of the 
state associated with band r is given by 



1 



\i+r 



sm 



ITTi 



N + 1 



(A6) 



and consequently the components of the unitary matrix which diagonalizes (^) are given 
by 



N + 1 



sm 



N +1 



(AT) 



The particle creation operators in the band and layer language are thus related by the 
equation 



TV 



ai 



(A8) 



Since we predominantly focus on intraband pairing the order parameter associated with each 
band r is thus defined through the anomalous Green's function components 



^r,r'(q) 



;ar-T(q)ar'i(-q) - ari(q)ar'T(-q)) 



N N 



N + 



-L 1 1 



sm 



i=i j=i 



ir-K 



N + 



Fi , (q) 6rr' 

N+1 '"^ 



(A9) 



We thus define the diagonal components of the order parameter in the band language 
by the equation 



N N 



A. 



A^ 



^ i=i j=i 



sm 



ITTl 

N +1 



sm 



N+1 



A 



(AlO) 



With the pairing restricted to intralayer pairing and pairing between nearest neighbour 
planes we can use the definition of Ai j to write the order parameter on a given band r in 
terms of the anomalous part of the Green's function 



N 



Ar(q) 



N + 



yE EMi(q,q') 



sm 



q' U=l 



N-1 

- 2 E ^J-lq-q) sin 

1=1 



N + 1 



N + 1 



sm 



ir TT 



sm 



^A^ + 1 
i + l)rT[ 



N + 1 



[F,,+i(q') + F,+i,(q')] (All) 



With a fully diagonalized Green's function the gap equation for the A^-layer system can now 
be readily obtained. Upon performing the sum over Matsubara frequencies in the usual way 
we get 



N 



A.(q) = - E E V^'^'^'(q,q')4^tanh 



q' r'=l 



2^.'(q') 



2/5^^'(q') 



(A12) 
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where the interaction in the band representation is defined as 

V^'-'iq, q') = ivTT [(l + I^^) V\\{ci, q') + (2 cos(^) cos(^) + 6;}/^) \4(q, q')] (A13) 



where the quasi-particle dispersions in the superconducting state are defined by the usual 
relation Ej. = \ ei + |Aj.|^. 



2. The Small t_L Limit 



For infinitesimal t±_ it is interesting to transform Eq. ( [A12| ) to the layer representation 
using Eq. (|A10|) . Then to zeroth order in t±_ the intralayer problems uncouple from the 
— 1 interlayer problems yielding the following two sets of gap equations 



tanh 



N 



(A14a) 



taiih 



|/5E(q'; 



N-1 



A...«(q) = ^I:V-X(q.q-) ' E (1 



+(N-1) 



A,,,+i(q') (A14b) 



Linearizing with respect to A and solving for the eigenvalues and eigenfunctions of these 
two equations gives the T^'s of the various pairing states available to the system and only 
one of these solutions will correspond to the maximal value of Tc. 

Let us first define two functions on the fermi-surface, ^/'||(q) and ilj±{q), and two numbers, 
fin and fl±, which satisfy the two equations 



^11 ^11 (q) 



dS 



Mi(q>q')^ii(q') 



(27r)2w^ 
r dS 



(Al5a) 
(Al5b) 



The integrations are over the degenerate bands of the Fermi surface. The pair amplitudes 
can be written in the form 



Ai,i(q) = Aji^ii(q) 
A,,+i(q) = A\^^{q) 



(A16) 



The gap equations (|A14| ) are transformed in the usual manner by separating the sum over 
quasi-momenta to separate integrals over energy and the Fermi surface. The energy integral 
can then be performed to obtain a BCS-like transition temperature and using Eqs. ( |A15|) 
we derive the linear matrix equations 



A^, 



In 



In 



keTj N + 



N-l 



ksTj N + 



(Al7a) 
(Al7b) 
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where 7 ~ 1.14 and uJc is the usual cutoff energy. We can now characterize all the available 
pairing states with functional form given by '?/'|| and ^l)^ along with their BCS transition 
temperatures. 



Eqs. (|A17|) each have one isotropic eigenvector given by = ^|j/±,o with Tc^ deter- 



mined by the equation 

i„(j^) = J- (AI8) 

which give the most stable candidate pairing states for the two mechanisms. The remaining 
eigenvectors correspond either to metastable states which are symmetric in the layer index 
such that Y.i ^||/± = with T^. given by 

or non-pairing states (T^ = 0) which are odd in the layer index. 

Transforming to the band picture using Eq. ( |A1(J| ) we find that the intralayer pairing 
states are even under the transformation r N + 1 ~ r whereas the interlayer pairing 
states are odd under this transformation. The two most stable candidate states thus give a 
superconducting order parameter of the form 

A,(q) = A||,o ^ii(q) - 2cos(^^^j A^.o ^±(q) (A20) 

The two parameters Ay and A_|_ o are both nonzero only when t±_ is finite and when ip^ 
and ijj± belong to the same irreducible representation. Since the two different pairing mech- 
anisms give solutions of different r dependence, the dominant type of pairing can be easily 
determined even when t± is finite. Mixing of solutions at finite t± and the effect on the 
transition temperature is discussed in Appendix P in the bi-layer context. 



3. Formulation In Terms of Bloch Waves 



For completeness we conclude by relating the above formulation to the usual treatment 
of layered materials. It would be natural to define a z-component of the quasi-momentum 
vector to be 



Qz 



N + 1 



(A21) 



Note, however, that there are linearly independent, non-degenerate eigenstates for 
< < vr and so taking —qz gives no new states. Thus interpreting q^ as a mo- 

mentum is rather unnatural. 

The usual procedure in the case of large A^, EE however, is to assume periodic boundary 
conditions. This means that we have hopping between layers 1 and A^ and that the system 
is translationally invariant along the c-axis. The eigenstates of this new Hamiltonian are 
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exp i 



2r7i 



with energies 



€r = 4 + 2t± COS( — 



(A22) 



(A23) 



Taking the hmit of A^ going to infinity one obtains the gap equation for a fully three dimen- 
sional system. Again assuming only intraband pairing the gap equation becomes 



A(q) = - ^V(q,q'),^i5lltanh 



2E(qO 



with the interaction 

^(q,q') 



A^ 



V||(q,q') +2cos(g^ -g^)\4(q,q'^ 



(A24) 



(A25) 



and the usual quasi-particle energies -E(q) = Ye(q)^ + |A(q)| . In the large A^ limit this 
system behaves identically to the one considered throughout this paper when only singlet 
pairing is considered. On the other hand, it is obvious that this formulation gives very 
different results in the small A^ limit. One might expect that in any real system the hopping 
between unit cells is different from that within a unit cell and so a more general formulation 
than either of these would be required. Such a formulation would, however, yield a contin- 
uous set of fermi levels within some small band and this result would be in contradiction to 
experimental observations. 



APPENDIX B: INTERBAND PAIRING 

To conclude this discussion we present a more careful treatment of the case for which 
N = 2 admitting the possibility of pairing of electrons on different sub-bands of the Fermi- 
surface. We will restrict our attention to states which are even under inversion, neglecting 
the possibility of a p-wave order parameter. This necessitates the consideration of a triplet 
interlayer pairing state. The normal state Hamiltonian has the form 

= ( 4. ) (Bl) 

where t± = |t±|e~'''''. The order parameter has the form 

A-( ^I + ^i Ai-^ELiAivA 

Ui + ^ELiAiV, M-A^ ] ^^^^ 



The (Tj are the three Pauli spin matrices. Upon diagonalizing ( |B1D we obtain 

- \tL\ 

^+\t±\ 

A+ Ai- 
Ai + 6A2 A 



A^f.^t._^';'^=) (B4) 



27 



REVTeX 3.0 



Levin Group Preprint, 1995 



We assume that the triplet component is describable by a single complex parameter and 
a real unit vector in spin space h so that A^* = A^nj. The Green's function can thus be 
diagonalized in its spin degrees of freedom. There are, therefore, four independent parame- 
ters which describe the superconducting state. The components of the order parameter in 
the band representation are 



(B5) 



Designating the bands by e-t = ^ =F \t±\ we solve for the anomalous parts of the Green's 
function as before. These may be linearized in A at and upon performing the sum over 
Matsubara frequencies we arrive at the following four gap equations 



A+ 


= Ajj + A 


A_ 


= Ajj - A 


Ai 


= A| 


A2 


= -Aids 



A, +A_ 



A, - A_ 



Ai 



A, 



q 
q 

EM|A 
q 

E^^^2 



tanh(i/3e+) tanh(i/5e_) 
A, Ai ^ + A_- ^ 



2e, 



tanhf i/3e+ I 
A+ Ai / _ A 



2e_ 
tanh(i/3e_) 



2e+ 2e 

tanh(i/5e+) + tanh(|/3e_) 
2e+ + 2e_ 

tanhQ/3e+^ + tanhQ/3e_^ 
2e+ + 2e_ 



(B6a) 
(B6b) 
(B6c) 
(B6d) 



If t_L = then there is a single transition temperature associated with the two intralayer 
pairing states and another transition temperature associated with the two interlayer states. 
Denote the larger one of these by Tc^. It can be shown that for non-zero t± the transitions 
described by Eqs. ( [B6c| ) and ( [B6cl| ) have a lower Tc than the intraband pairing state transition 
in Eqs. ( P6a|) and (|B6b|) where for small t± the new transition of the interband pairing stateil 
with higher is given in terms of T^, by 



In 



-0.2123 



The transition temperature for the intraband pairing states is given by 



In 



T, 



c,0, 



In 



t±R 



(B7) 



(B8) 



where the overlap between the two pair wavefunctions is 

dS 



R 



{2'kYvf 



(B9) 



We thus see that in the presence of interlayer tunneling the favoured pairing state is always 
an intraband pairing state. 
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TABLES 







J. 

H 












i 


U 


U 


U.ioUo 


i ( COS On^ -\- COR Ooil 


-1.0 


-0.50 


-0.5951 


cos Qx COS 


0.2 


0.15 


0.1636 


1 (cos 2qx + cos 2(/y) 





-0.05 


-0.0519 


^ (cos 2qx cos + cos Qx cos 2gy) 








-0.1117 


cos 2qx cos 2qij, 








0.0510 


i (cos qx - cos 0^2/) 





0.04 
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0.08 


TABLE I. Tight binding basis functions and hopping parameters (in eV) used 


in numerical 



calculations. 
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FIGURES 

FIG. 1. Phase diagram showing order parameter symmetry in LSCO: (a) Near van Hove points; 
(b) Away from van Hove points. The two axes represent the wave vectors Qx, Qy at which the 
pairing interaction has a maximum. 

FIG. 2. Phase diagram for one layer model of YBCO: (a) Near van Hove points; (b) Away from 
van Hove points. The dotted line separates Aig representations with and without a, 6-axis tt phase 
shifts. 

FIG. 3. Phase diagram for one layer model of BSCCO: (a) Near van Hove points; (b) Away from 
van Hove points. A dotted line separates states with four and eight lobes within a representation. 

FIG. 4. Phase diagram for stability of vr phase shifts on the a relative to the b axis, in the 
order parameter (below solid line) Model is for YBCO with an attractive interlayer and repulsive 
intralayer interaction, both peaked at Qaf- Hole doping fraction appears on the ordinate and 
relative strength of the two interactions defines the abscissa. The shaded regions denote states 
where the order parameter on one (light) or both (dark) bands is nodeless. 

FIG. 5. Evolution of order parameter solutions as interlayer correlations are increased. The 

states correspond to the solutions obtained in the previous figure for x = 0.25 and the indi- 
cated A_l/A||. Observe that proximity to the van Hove singularity results in considerable a, 6-axis 
anisotropy despite the very small orthorhombicity. 

FIG. 6. Schematic illustration of Josephson coupling between the anti-bonding band in YBCO 
(left) and a Pb counterelectrode (right). Overlap of the single particle wave functions is nonzero 
unless the two wavefunctions are perfectly aligned in the directions parallel to the plane of the 
interface. 

FIG. 7. Possible scenarios for order parameter behavior across twin boundaries showing tt phase 
shifted (c?, d)-type (a) or (s,— s)-type (b and c) solutions. The scenarios depicted in (a) and (c) 
give TT junction behaviour in a, 6-axis corner junctions while case (b) will lead to cancellation of 
the TT phase shift after averaging over twin domains. 

FIG. 8. Amplitude of pure intralayer and interlayer pairing solutions in the band representation 
for = 4 (a and b respectively ) and = 7 (c and d respectively) cases. The ideal solutions 
{t± = 0) are denoted by the dotted lines and the actual solutions are calculated for moderately 
large t±. 

FIG. 9. Solution for a three band problem illustrating co-existence of solutions of different 
symmetries. The inset plots the associated band structure. 
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FIG. 10. Pairing Density of States (PDOS) of the solution with highest Tc in each irreducible 
representation for LSCO. Panel (a) shows the absolute PDOS in arbitrary units and panel (b) 
illustrates the PDOS normalized by the single particle density of states. The isotropic s-state (solid 
line) was calculated for an attractive pair interaction peaked at (tt, vr) while the other solutions 
were calculated for repulsive interactions, with Q chosen so as to maximize Tc- 
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